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ABSTRACT 

Multiscale phenomena are ubiquitous in nature as well in laboratories. A broad 
range of interacting space and time scales determines the dynamics of many systems 
which are inherently multiscale. In most research disciplines multiscale phenomena 
are not only prominent, but also they have often played the dominant role. In the 
solar wind - magnetosphere interaction, multiscale features coexist along with the 
global or coherent features. Underlying these phenomena are the mathematical and 
theoretical approaches such as phase transitions, turbulence, self-organization, frac- 
tional kinetics, percolation, etc. The fractional kinetic equations provide a suitable 
mathematical framework for multiscale behavior. In the fractional kinetic equations 
the multiscale nature is described through fractional derivatives and the solutions of 
these equations yield non-convergent moments, showing strong multiscale behavior. 
Using a Levy-flights approach, we analyze the data of the magnetosphere and the 
solar wind. Based on this analysis we propose a model of the multiscale features and 
compare it with the solutions of diffusion type equations. The equation with frac- 
tional spatial derivative shows strong multiscale behavior with divergent moments. 
On the other hand the equation with space dependent diffusion coefficients yield 
convergent moments, indicating Gaussian type solutions and absence of long tails 
typically associated with multiscale behavior. 
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1. INTRODUCTION 

Many laboratory and natural systems have a broad range of interacting space and 
time scales and their dynamics exhibit multiscale features. The time-series data of 
such systems have some common characteristics: the physical variables reveal the 
properties often referred to as multiscale phenomena, heavy tails, strong intermit- 
tency, etc. Examples of such systems include dynamical systems with chaos, stock 
market indices, earthquakes, space and laboratory plasmas, atmospheric and hydro- 
logical systems, etc. A variety of characteristic quantities, such as Hurst exponents, 
multifractal spectra. Levy flights, etc. have been computed from the time series 
data using different approaches and techniques. Typically, the application of a spe- 
cific method has been largely motivated by the current understanding of the relevant 
processes in a particular system, and their relationship to the available data. 

From the dynamical systems point of view, the multiscale properties suggest frac- 
tal or multifractal structure of the dynamical processes and numerous simulations 
support this idea [1]. The advantage of the analysis of data generated by chaotic dy- 
namics is a possibility to link the system evolution to the equations of the processes 
and to their origins, e. g., physical, chemical or biological processes. For example, the 
Levy flights can be identified as dynamical processes of trapping into potential wells 
or ballistic propagation. Observing and analyzing different realizations of stochas- 
tic dynamical processes, one can classify them and apply the available techniques to 
the phenomena of which the details of the processes are not sufficiently known. A 
technique similar to the one used in the analysis of Levy flights has been developed 
recently and has been applied to the plasma density fluctuations in the tokamak de- 
vice [2]. Many studies of the solar wind - magnetosphere coupling have used the 
time series data to obtain the Hurst exponent characterizing the multiscale nature 
[3,4,5]. These studies have motivated the modeling of multiscale behavior in terms of 
well-known equations such as Fokker-Planck [6, 7] and Ginzburg-Landau equations 
[8]. 
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It is however very desirable to apply different techniques to the same data set 
to obtain a consistent characterization of the underlying processes, and thus provide 
a strong basis for developing models. This point is underscored by the fact that 
the knowledge of the Hurst exponent [9, 10] does not provide sufficient information 
about the physical processes underlying the data. Indeed, the expression, < [x{ti) — 
x{t2)Y >~ 1^2 — which defines the Hurst exponent H, in fact, does not provide 

a clear and unique information about the random process x{t). This is an important 
issue and three pertinent points are elaborated in the following. 

First, the data typically characterize some specific physical observable x{t) and 
the Hurst exponent is strictly related to that observable. Other physical variables will 
have different exponents and the connection between these exponents is not known 
in the absence of more information about the real processes. Practically, it means 
that the common characterization as "sub diffusion" for if < 1/2 and "superdiffu- 
sion" for if > 1/2 have no real meaning unless the physical processes underlying the 
diffusion or, more generally, the kinetic processes, are specified. Second, the Hurst 
exponent is related to the second moment only, and does not reflect the behavior of 
other moments. For example, for Gaussian or Poissonian processes the second mo- 
ment is sufficient to characterize the system as all higher moments can be expressed 
in terms of it. Some Gaussian type processes can have H < 1/2, and the higher 
moments can be evaluated from the data that yield the Hurst exponent. However, 
this feature can be revealed only by computing the higher moments. For the Levy 
and other multiscale processes the second and higher moments are infinite. Such 
features are not captured by the Hurst exponent even if if > 1/2 since, as mentioned 
above, the process may be super diffusive, but it can not lead to a definite identifi- 
cation. Third, the Hurst exponent represents the time dependence of the variance 
but the time characteristic of the process, typically, is not sufficient, as is clear in 
the case of chaotic dynamics [1]. Indeed, the physical process can be characterized 
by a probability distribution function F{x,t) and the moments of x are defined as 
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< x'^{t) >— J x'^{t)F{x{t),t)dx{t). In order to distinguish, for example, a multiscale 
process from a Gaussian one, more subtle features from the data than just the Hurst 
exponent arc needed. It should be noted that even for a purely diffusive process for 
which F{x, t) satisfies a kinetic equation with derivatives with respect to x and t, not 
only the exponent H is important for the kinetics but also the deviation 5x during 
some specified interval 5t is required. 

The important role of kinetics was effectively demonstrated in [zaslavskyOO, 2] 
where a new technique of data analysis was applied to the data of plasma density 
fluctuations in a tokamak. Assuming the presence of a Levy-type multiscale process 
characterized by "flights" , a special procedure was used to separate the flight events 
without a characteristic scale from the "noisy" events with a characteristic small 
scale. The notion of "flight" is however non-trivial and it should be defined carefully 
to avoid ambiguity. It should be pointed out that recent studies of the data of space 
plasmas [6, 7] did not define the Levy flights in their analysis. This issue deserves 
more careful attention [1]) and will be discussed in more detail in the following. 

A simple way of introducing the notion of a Levy-type flight, or in short, flight, 
is to consider a function x(t) or its integral(s) and obtain the long intervals St where 
the changes in x{t) or its integral (s) are monotonic up to some level of accuracy. The 
change Sx in x{t) is then a measure of the length of a flight with the distribution 
function F{5x,5t), within the specified level of accuracy. It is very important to 
note that F{6x, 5t) is the distribution of flights if they exist! Further, all moments 

< x'^{t) >, obtained from the data, are always finite, i.e., F{Sx,St) is truncated and 
defined up to Sx„iax, ^tmax- Then the extrapolation of F{5x,5t) should be carried 
out with utmost care. For example, in the case of F{Sx, St) ~ l/{Sx)'', the moments 

< {6xy > with q > k — 1 will be infinite, but the observational data can yield the 
exponent k within the interval < 5x < Sxmax- 

In this paper a specific form of Levy fiight analysis [2] to the data of the solar wind- 
magnetosphere coupling. A correlated dataset [11] of the magnetospheric dynamics 
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driven by the solar wind is used in this study. In this dataset corresponding to the 
last peak of the solar cycle, the solar wind (SW) is represented by the flow-induced 
electric field and the magnctospheric response is represented by the auroral electrojet 
index AL. This data set is used to study the multiscale properties by applying a 
technique based on an a priori assumption of the existence of Levy-like flights. This 
analysis yields a new characterization of the magnctospheric dynamics. 

Multiscale behavior has been modeled using the so-called fractional generalization 
of the Fokker-Planck equation [1]. In order to determine the multiscale properties of 
such equations, numerical solutions of specific cases of such equations are studied. 
In the case of scale dependent diffusion coefficients it is found that the solutions of 
the equations yield convergent moments, and thus are not suitable for modeling the 
multiscale properties. In the case when the equations have fractional derivatives, viz. 
they represent fractional kinetics, and the numerical solutions do not yield convergent 
moments. It should be noted that the diffusion equation, even with a fractional time 
derivative representing the so-called fractal Brownian motion, has solutions in the 
form of a Gaussian function G{x/t^), for which all the moments are finite for any 
exponent a. In the case of a kinetic equation with a fractional derivative with respect 
to coordinate or momentum one can get the same exponent but only few (or even 
a single) first moments are finite, and all higher moments are infinite for any t. In 
the case when the kinetic process F{x, t) is modeled using the Fokker-Planck type 
equation with a scale dependent diffusion coefficients, as in [7], we show that all 
moments are finite. This leads to the conclusion that a correct exponent for the 
second moment alone is not sufficient to construct a model, and it is essential that 
the model yields the appropriate multiscale behavior. 

The structure of the tails of the distributions of different physical observables 
becomes crucial (see also [12]), and as mentioned above, for data of finite length, 
defining the tail is a fairly sensitive task and needs confirmation from independent 
measurements. In particular, a specific demonstration is provided on how the accu- 
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racy of the data for the tail can be improved (see Fig. 4) if the relevant processes 
exhibit flights. 

The structure of the paper is as follows. Section 2 presents a brief overview of the 
solar wind - magnetosphere interaction, the database and the recent studies relevant 
to this paper. In Sec. 3 the Levy flight technique used to analyze the data [2] is de- 
scribed and also presents the corresponding results. The analysis indicate a power-law 
type behavior of the distribution functions in time and in the magnitudes of the fluc- 
tuations. The multiscale character of the observational data and some mathematical 
models for their description are discussed in Sec. 4, followed by a detailed analysis of 
the fractional kinetic equation in Sec. 5. The reconstruction of the kinetic equations 
from the data is discussed in Sec. 6. Numerical studies of the fractional kinetic equa- 
tion is described in Sec. 7, and the conclusions of the paper are summarized in Sec. 8. 

2. COUPLED SOLAR WIND - MAGNETOSPHERE SYSTEM 

The Earth's magnetosphere is a huge cavity, mostly shielded from the flow of 
charged particles coming from the Sun, the solar wind plasma, by its magnetic field. 
The solar wind plasma can enter this magnetic shield through magnetic reconnection. 
The main source of plasma for the magnetosphere is the ionosphere, a relatively thin 
and dense plasma shell, which separates the magnetosphere from the Earth's thermo- 
sphere and atmosphere. The solar wind flow is not steady, because of the active solar 
corona, e. g., solar flares, coronal mass ejections and other bursty processes, as well 
as the inherent turbulent nature of the solar wind flow towards the Earth. Similarly, 
the magnetospheric activity is determined both by its interactions with the solar wind 
and by the complex internal activations on many scales occurring in the practically 
coUisionless plasma of the magnetosphere. Attempts to classify the variety of ob- 
served phenomena have resulted in descriptive classifications such as magnetospheric 
storms, substorms, pseudo-breakups, convection bays, saw-tooth events, etc. Yet, the 
description remain essentially incomplete. Moreover, with some exceptions, e.g., the 
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storm-substorm relationship, they do not reflect important connections and similar- 
ities between processes occurring on different spatial and temporal scales, which are 
revealed both in the solar wind [13, 14] and in the magnetosphere [15, 16]. 

The recent studies of the multiscale behavior in the coupled solar wind - magne- 
tosphere system have been based on two approaches. In the flrst, the observational 
data from different sources, such as spacecraft-borne and ground-based instruments, 
have been used to analyze the scaling properties. In the second approach, mathemat- 
ical models such as Fokkcr-Planck type equations have been used to interpret these 
processes. In the latter equations which are known to exhibit multiscale properties 
are used to study the conditions under which the scaling laws and exponents derived 
from data can be obtained. However it should be noted that the ability of a partic- 
ular type of equation to flt the data does not exclude other types of equations from 
reproducing the desired scaling, etc. Such studies provide valuable results that may 
lead to the origins of the multiscale behavior. These two approaches, viz. the studies 
using data and the modeling using equations, are discussed in the following sections. 

The magnetospheric activity is characterized by a set of parameters, including 
both the basic physical ones, such as the magnetic and electric flelds, plasma density 
etc., and the specific geophysical parameters or indices, e. g., Kp, Dst, AL, AU, AE, 
etc. [17]. A particularly relevant set of indices is the family of Auroral Electrojet in- 
dices AE, AU and AL, obtained from a number (usually greater than 10) of stations 
distributed in local time in the latitude region that is typical of the northern hemi- 
spheric auroral zone. For each of the stations the north-south magnetic perturbation 
H is recorded as a function of the universal time. A superposition of these data from 
all the stations yields a lower bound or maximum negative excursion of the horizontal 
component of the magnetic field variations, and is called the AL index. An upper 
bound, defined similarly, is called the AU index and the difference between these 
two indices, AU-AL, is called the AE index. Auroral indices are particularly relevant 
for substorm studies and have been used extensively for more general nonlinear time 
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series analysis and forecasting [18, 19] and for multiscale studies [16]. 

The scaling properties inherent in a data set can be represented by many character- 
istic parameters. One of the most popular is the Hurst exponent H [9, 10] measuring 
the "roughness" of multifractal data and allowing one to distinguish between con- 
ventional diffusion (H=0.5), subdiffusion {H < 0.5) and superdiffusion {H > 0.5), 
as discussed earlier in Section 1. However, the presently available results for the 
Hurst exponent, computed largely from the AE index are rather controversial, can 
not distinguish between sub- and super-diffusion regimes. In particular, Takalo and 
Timmonen [3] found H=0.5 based on the AE index studies. Uritsky and Pudovkin 
[20] concluded that H is largely less than 0.5 (H=0.38-0.59, depending on the sub- 
storm phase). In contrast. Price and Newmann [21] found that H is mostly greater 
than 0.5 (H=0.44-0.97, depending on the time scale range). In a series of studies 
involving extensive data processing, Hnat et al. [4, 5, 22] used different geophysical 
variables, including the solar wind parameter e, all three auroral indices, and the 
square of the interplanetary magnetic field S^, to obtain the inherent scahng. These 
studies showed most of the Hurst exponents to be subdiffusive [7], with H ranging 
from 0.28 to 0.52. This is the most reliable result and it confirms the fractal nature 
of the solar wind. Overall, the observed data suggest subdiffusive behavior, which 
is drastically different from the superdiffusive behavior demonstrated by the widely 
used anomalous diffusion models utihzing the concept of Levy flights [1, 2, 23]. 

Attempts at resolving this ambiguity [24, 25] (see also [4, 5]) point out that only 
some of the parameters reveal the actual scaling corresponding to monofractals. A 
good example is for the solar wind parameter [4]. In other cases either the gen- 
eralized structure function (GSF) (depending on its order m) or the probability dis- 
tribution function (PDF), cannot be rescaled. Moreover, typically the effective H, 
inferred from the GSF of the order m=l, is more than 0.5, suggesting a superdiffu- 
sion, while H inferred for m > 1 is equal to or less than 0.5, suggesting a subdiffusion. 
Another reason mentioned by Watkins et al. [25] is that almost any simulation of 
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fractional Levy motion is effectively one of truncated Levy motion, which may have 
the C parameter of the GSF superdiffusive corresponding to the case for m=l and 
diffusive (H=0.5) for m=2. 

Thus, the GSF analysis shows that in many cases the solar wind and magneto- 
spheric data do not show ideal scale invariance, that is, they do not correspond to a 
monofractal. And as a result, the simplest conclusions in terms of the Hurst expo- 
nents are too ambiguous and are not clear enough even to distinguish between sub- 
and superdifFusion behavior. On the other hand, Hnat et al. [5] revealed a collapse of 
the rescaled PDFs for the solar wind parameter e and the AE index, suggesting that 
some form of simple scale invariance may exist and may even allow one to distinguish 
between these scaling properties for the solar wind and the magnetosphere. Accord- 
ing to Hnat et al. [22], this does not exclude the possibihty that these scaling effects 
can be explained in terms of classical Fokker-Planck diffusion with scale-dependent 
diffusion coefficients. Finally, a fairly general analysis of the data [7] increases the 
level of ambiguity in the interpretation of the data. 

The pitfalls of the data analysis in the case of multiscale processes, discussed in 
Section 1, leads to a couple of comments on the results from some recent studies [7, 
8, 22]. First, in these papers the dispersion of the data in the tails is too high to draw 
a clear conclusion on whether a power-law tail exists. Second, the claim of a possible 
multifractality is not specified and, in fact, such behavior has two origins. One is that 
the data can have a few different scales (as it will be seen from our analysis in the 
following) in different intervals of t. This may not be a surprise since different physical 
processes may be responsible for different intervals of the magnetospheric data. The 
second source of multifractality is related to the so-called log-periodicity [1], i.e. a 
small modulation of the scahng law in logt scale. This effect, while not discussed in 
[7, 8, 22], can be seen in Figs 1, 3, 4 of [7]. The existence of log-periodicity requires 
very different data analysis [1] in which the behavior of the higher moments becomes 
crucial. 
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The alternative approach to the data analysis and their interpretation, presented 
in the next section, is aimed at enhancing the information extracted from the obser- 
vational data, and at improving the accuracy of its interpretation. It is worthwhile to 
add that this approach and technique have been developed building on the experience 
in the study of different models of dynamical chaos. The study of these models with 
the new technique has led to the improvements in its accuracy and has enabled one 
to verify some a priori assumptions on the character of processes with flights. 

3. ANALYSIS USING LEVY FLIGHT TECHNIQUE 

Let y be a physical variable of interest, for example, the density of ions, or a 
component of the magnetic field, or the plasma current in a specified direction, and 
let Ui be its value at discrete time instant ti. The time instants {ti, . . . , t„} are taken 
to be equidistant. It is also assumed that all data {yi} are collected at the same 
location, i.e., the sequence {yi} is a pure time-series at a specified point in space. 
Alternatively, the time series can represent the global behavior, e.g., the indices in 
the case of the magnetosphere. 

The integrated value 

S(t) ^ f dty{t), (1) 
Jo 

or, for a discrete time series {ti}, 

Sn^SiQ^f^Vr (2) 

is used in the analysis of flights. In the case when the law of large numbers is valid, 
the fluctuations of Sn with respect to its mean value (Sn) satisfy the condition 

{\5S{tj,)\') = - (Sj,)]') ~ TV , TV » 1, (3) 

as in the case of Gaussian processes. In the anomalous case with power-law tails of 
the distribution of the values 5'„, the relation (3) does not hold since the second or 
higher moments of the probability distribution function of A^, P{AS), are inflnite 
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with 



ASn = Sn- n{yn - yo) 



(4) 



The corresponding processes are of Levy type [26, 27]. The new variable ASn in- 
troduced in (4) differs from the integrated variable -S'^ by the subtraction of a hnear 
trend. This replacement is done for convenience in working with the data, and it does 
not influence any physical interpretation of the results. 

The solar wind induced electric field and auroral electrojet index [11] are used 
to compute A5'sw,n and A5'AL,n and the results are shown in Fig. 1. The results 
show large and small scale wiggles. A few zooms in Fig. l(c,d) show a kind of self- 
similarity of the behavior of ASn, ic, an indication of the presence of multiscales. 
Further analysis of the data can be done using a kind of pattern recognition scheme. 
Consider a set of connected values in Fig. 1 as a curve with many different segments 
and each segment is defined as an interval of a monotonic behavior of the curve ASn- 
Monotonicity of the curve ASn within an interval Aua — {n^^\n^'^^) will be defined up 
to some interval a which is the interval of coarse-graining. This means that possible 
values of Aua form the set 



The monotonic piece of the curve ASn within an interval Aua will be called a 
"flight" and Aua is the duration of the flight. All flights are defined up to a smoothing 
interval a. The flights can be deflned in both directions, i.e., as monotonic increases, 
or decreases of ASn within an interval of monotonicity. 

The set of data of length ~ 10^ is used to compute the probabihty distribution 
function Pgwijia) and PAL(^^a), i-e- the distribution of the time duration of flights 
for the solar wind and auroral lower index, respectively. The corresponding results 
are shown in Fig. 2 for three different values of a = 2,4 and 6. It is evident that 
the resulting curves do not differ strongly. The results for n2 (a = 2) are chosen for 
detailed analysis, without a loss of generality. 



Aua e {Ana} = {o, 2a, 3a, . . .} 



(5) 
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Figure 1: Coarse-grained data time series for the deviation ASsw,n (Solar Wind) and 
^SAL,n (Auroral Lower index) as a function of discrete time n: (a) real data (solid 
line) for ASsw,n and refined data with extrapolation of the data values for data gaps 
(dot-dashed line); (b) real data for ASAL,n', (c) and (d) - zoom of the corresponding 
curves in (a) and (b). 
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Figure 2: Plot of the distributions P{An) of the time-length An of "flights" for the 
solar wind (SW) and the magnetosphere (AL), with averaging over a — 2 (a), and 
a — 6 (h) points along n. 

Both the distributions Psw(^) ^-nd PAhin) have three characteristic segments de- 
fined by intervals of n: zone 1 (0 - 10), zone 2 (10 - 300), zone 3 (> 300). In the first 
and third intervals the distributions can be written as 

P(n) ~ const/n^+^ (6) 

where the exponent /3i depends on the type of data and on the interval. The computed 
values of the exponents are 

P^aI ~ P^w ~ 0.9 - 1.0 
/3S ~ 1-3 - 1.8 

P§i, ^ 1.5 - 2.2 (7) 

The big error bar for zone 3 is due to an ill-defined width of the transition zone 2. 
The slope corresponding to AL appears to be smaller than that for SW in zone 3: 

Pal < Psw^ (8) 
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and the accuracy of the computed values is higher, the larger is the value of a for the 
selected intervals of coarse-graining. 

The transition zone 2 looks smoother for the solar wind data and one can assume 
that this zone corresponds to some resonant type process whose origin is in the solar 
wind. The response of the magnetosphere to the corresponding interval of flights 
duration makes the transition zone sharper. 

Similarity of the behaviors of SW and AL data in zone 1 makes it possible to 
predict short flight evolution in the magnetosphere. A similar possibility to predict, 
though with a different probability, can be applied to zone 3. 

There is a possibility to extract additional information from the data in Fig. 1. 
It is the change in the magnitude of Sn during a flight An^. Let s be such a change 
and Psw(s), -Pal(s) the probability distribution function for such changes. One can 
say that s is a change of the corresponding magnitude for an intermittent process 
during the period of monotonicity. Following the definitions (1), (2), and (4), during 
the monotonic flight between time instants ti, t2 the value of s is 



corresponding to the so-called parabolic flights with a constant acceleration. The 
corresponding distributions are given in Fig. 3 for two cases of coarse-graining: a = 2 
and a — 4, which do not show a difference. The results in Fig. 3 provide interesting 
observations. There are no 3 zones as in Fig. 2. The level of changes of s in SW is 
larger by almost an order of magnitude than in AL as it is seen from Fig. 1. The 
same difference is seen for the distribution functions from Fig. 3 within interval of 
s e (10 - 300). After s > 1000 the PDF has the form 



Si2 = const. (^2 - ^i) , 



(9) 



corresponding to the so-called ballistic flights, or 



Si2 = const. (^2 — ^l)^ 



(10) 



P(|s|) ~ const/|s 



l+ai 
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with ai given by 

ttsw ~ 0.3 — 0.6, 

ttAL ~ 0.7 - 1. (12) 

The most interesting part of the distributions in Fig. 3 is that for s ~ 1000, the 
difference between PALd-sj) and Psw(|s|) is almost negligible for s e (10^ - 10^). This 
indicates that in this range of s the magnetosphere and the solar wind have similar 
characteristics. The studies of burst life time distributions [28] are closely related to 
these studies. 

The above values of a correspond to the super diffusive case. This is in contrast to 
the conclusion drawn from the values of the Hurst exponent, which has values < 0.5 
and thus corresponds to the subdiffusive case. As mentioned in the introduction, the 
explanation of the Hurst exponent needs more specific details of the physical pro- 
cesses and the corresponding theory (Sees. 4, 5; see also discussion in [4, 5, 25]). It 
should be noted that only some of the parameters reveal the actual scaling, corre- 
sponding to monofractals. A good example is the case of in the solar wind [4]. 
In other cases either the generalized structure function (depending on its order m) 
or the PDF cannot be rescaled, or both. Moreover, typically (say, for AE) the effec- 
tive H, inferred from the m = 1 GSF, is more than 0.5, suggesting a super diffusion, 
while H inferred for m > 1 is equal to or less than 0.5, suggesting a subdiffusion. 
It may be noted that a truncated Levy motion may also have the ( parameter of 
the GSF corresponding to superdiffusive for m=l and diffusive (ife//ectTOe=0.5) for 
171 — 2. The inconclusive nature of these studies was discussed in the previous section. 



4. MULTISCALE FEATURES IN OBSERVATIONAL DATA AND THEIR 
MATHEMATICAL MODELS 

Many systems such as the coupled solar wind - magnetosphere system described 
in section 2 exhibit multiscale behavior in the form of power-law distribution of scale 
sizes or long tails in the probability distribution functions. The mathematical frame- 
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Figure 3: Plot of the distributions -P(s) of the change s in the magnitude Sn per 
flight for the solar wind (SW) and the magnetosphere (AL). 
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work for modeling these phenomena, based on the diffusion-type equations, are de- 
scribed in this section. 

For a system that exhibits self-similarity, the probability distribution function of 
a variable y{t) can be expressed as 

F{y, t)dy = F{y/ndy/t^ = F{Od^ , = y/n. (13) 

The exponent u is to be determined from the dependence of the first moment on t as 

{\y\) = J \y\F{y, t)dy = const.i^ (14) 

with the symmetric F{y,t) and the normalization condition 

jF{y,t)dy^l. (15) 

For example, if y is a coordinate that is subject to a diffusive type random walk 
process, then 

H^2v (16) 
and is called the transport exponent since for the second moment 

(7/2) = const.t^ . (17) 

The expression (17) needs further comments. The regular, or normal, diffusion 
has 11 = 1. Different values of /i, including = 1, can result from the Fokker-Planck 
type equation 

dF{y,t) d ,dF(y,t) 

depending on the function T>[y) that describes nonuniformity of the diffusion coeffi- 
cient [29]. The main feature of the diffusion process is that F{^) decays exponentially 
as ^ — > ±oo and, as a result, all moments are finite: 

{y"") < oo , 0<t <oo (19) 

for all n. 
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Another situation, when Eq. (19) is vahd is the so-called case of sub-diffusion [30] 
described by 

d'^F(y,t) d ,dF(y,t) , , 

where 

< /3 < 1. (21) 

In Eq. (20) the time variation is described in terms of a fractional derivative of order 
/3, and this is an example of a fractional kinetic equation. This and other forms of 
fractional kinetic equations will be discussed in more detail in the following. For 
T>{^y) = Vq = const., Eq. (20) gives 

(y^) = const.i^, (22) 

which is slower than normal diffusion as a result of Eq. (21). In the both cases of 
Eqs. (18) and (20) F[y,t) have a characteristic width 

AC = Ay/t'^ (23) 

such that for ^ ^ A^ the values of F(^) are exponentially small, which imphes finite 
moments {y'^) at any finite t. 

The situation is very different when the probability distribution function F(|/, t) 
behaves as 

cit) 

Fiy,t)^^, <^>0, t^oo, y^oo, (24) 

y 

where c{t) is a function of t alone. Then exponent 5 can be obtained as the slope in 

a In F{y,t) vs. Iny plot. An important part of the data-derived models is a plot of a 
histogram of this dependence. The elementary connection between variations 

Alny = {l/\y\)Ay (25) 

for y ^ oo shows that a collection of the data of In F{y, t) into equal bins along the 
Iny axis would provide less dispersion than the collection of the data into equal bins 
along the y axis. This statement is demonstrated in Fig. 4, where the log of the 
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probabilities P{An) ior a — 2 is plotted as a function of InAn with bins of equal size 
in An. 

Another, equally important comment is related to the truncation of F{y, t) for 
fairly large y, i.e., one should consider the data for F{y,t) and truncate part of the 
data for \y\ > ymax with t/max that is not well defined. Let the experimental or 
simulation (called "raw") data for a fixed t be within the interval y e (yo,yi) and 
\yi\ ^ \yo\- The power law distribution (24) imposes large fluctuations of the data 
for large \y\ with a finite (though not exponentially small) probability. Independent 
of how big t is and how much larger yi is compared to y, there will always be a lack 
of statistics close to yi, and this is why some data within an interval (ymax, |yi|) has 
to be excluded. 

There is a complementary result for the truncation of the raw data in all the 
moments 

(|2/r)tr<oo, (Vm), m>0 (26) 

are finite and this permits the computation of the self-similarity and even the multi- 
fractality parameter C{^) as 

(|yr)tr~^^^"^^ (Vm), m>0 (27) 

for the function F{y,t) if it exists. In the monofractal case 

C(m) = mii/2, (28) 

as compared to Eq. (17). 

In reality, the mult i-f r act ality can be due to the so-called log-periodicity [1], as 
discussed earlier. In particular, this means that there is some ambiguity in obtaining 
C(m), /I, S, etc., which depend on ^max- Thus for any exponent obtained from data, 
the interval for which it is computed and the quality of the statistics should be spec- 
ified. Comparison of the exponents from different data sets must be within similar 
intervals of the definition of the exponents and for similar statistics. 
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Figure 4: Log-Log plot of the distribution function P{An) with bins of equal sizes 
(equal to 2) in An, as compared to Fig. 2 where the bin sizes are equal in log An. 
Dots are for the AL data and circles are for the SW refined data. The plot indicates 
a strong dispersion of the distribution, in comparison with the plot in Fig. 2. 
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5. FRACTIONAL PHENOMENOLOGY AND FRACTIONAL 
KINETICS 

The crucial feature of the probabihty distribution function F{y^ t) is the existence 
or non-existence of the stationary distribution F{y) for t oo. In many physical 
problems there is no stationary distribution and the phase space of the problem has 
at least one unbounded variable. The plasma density fluctuations in laboratory fusion 
experiments (tokamaks) have been modeled in terms of a fractional kinetic equation 
[31]. Another example is the multiscale nature of the magnetosphere where the corre- 
lated data of the solar wind - magnetosphere system shows clear non-Gaussian PDF 
[32] . Zaslavsky [23] proposed a phenomenological approach to the nonstationary case 
of the evolution of t) in the presence of multiscale processes. 

The basic approach is to write down a balance equation 



5tF{y,t)^SyF{y,t\ (29) 

where the bar implies averaging over all admissible paths. This equation can be 
written as 

StF{y,t) = A^F{y,t) (30) 

where Af is a generalized difference operator for a time shift by At parameterized 
by P (see below). Similarly 

5yF{y, t)^Y: ^^['4(Ay, y)F{y, t)] (31) 

Ay 

where is a generalized difference operator for a y- variable shift by Ay along a 
path A{Ay, y) in the phase space, parameterized by a, and the line in (29) implies 
averaging over all admissible paths. Definitions of Af and A^ are given in [1] and 
their main properties are 

Af ~ {Atr^ , {At ^ 0) 



^y-i\^y\r^, (|Ay|-0) (32) 



22 



where the fractional derivatives of the orders (5 and a, respectively, are introduced. 
These derivatives work in different ways since there are different definitions of the 
variables t G (0, oo) and y G (—00, 00). 

To simplify the problem we consider t) to be symmetric with respect to y 
and A — A{/^y). Then the balance equation can be written in the form [zaslavsky02, 
zaslavskyOO; 1, 23]: 

where 

2^a/3 = ^-^A{Ay) , ^(0) = const. (34) 

Eq. (33) is a fractional kinetic equation (FKE) [1, 23] since (a,/?) can be fractional. 
In the case of f3 = 1, a = 2, Eq. (33) is a regular diffusion equation with the diffusion 
coefficient Va/s- The important part of the derivation of (33) is the existence of finite 
T>ai3 in the limit 

At ^ , Ay ^ , Vap = const. , (35) 

known, for the case /3 = 1, a = 2, as the Kolmogorov condition [1]. The existence of 
the limit (35) is a nontrivial fact, and eventually it defines the variable y for which 
the FKE can be written. 

An important property of Eq. (33) and the definition (34) is the rescaling invari- 
ance. Consider the following rescaling 

Ay ^ A^Ay , At ^ XtAt (36) 

that gives 

^ {X^/X?)V^p . (37) 
The FKE Eq. (33) is scale-invariant under the transforms (36) and (37) if 

p/a = lnXy/hiXt . (38) 

On the other hand, multiplying Eq. (33) by and integrating both sides over y, 
yields 

(|yr)tr = const, , (39) 
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i.e. the transport exponent 



II — 2P/a, 



(40) 



or using (38) 



II — 2ln Xy/ In . 



(41) 



For P = 1, < a < 2, the FKE describes the Levy process. The moments {\y\^) 
are finite ii p < a [30] and Eq. (39) can be considered only for the truncated F{y, t) 
as was commented on earher. 

In simple terms, if the data of a physical system exhibit the scaling property given 
by Eq. (36), then the FKE, Eq. (33), can be considered as a model of the processes 
the data represent. 

6. FKE FROM THE DATA 

The modeling of a physical system, whose data yields the scaling properties de- 
scribed in Sec. 2, in terms of a fractional kinetic equation is discussed in this section. 
There exists an indefiniteness in the reconstruction of the parameters {a, (3) that de- 
fine the FKE. Starting with the rescaling parameters {Xy, At) obtained from data, only 
the ratio (3 /a, but not the independent values of {a, (5) can be obtained. However, this 
indefiniteness does not infiuence the transport exponent ii. Another indefiniteness of 
the data analysis arises from the Kolmogorov condition (35) since this information 
can not be obtained from the the time series data alone, and other approaches need 
to be employed. 

In Sec. 2 the data of the solar wind - magnetosphere system were analyzed to 
obtain the distributions P{s) and pit). Since s is the accumulated variations in the 
variable, for example, the magnetic field in the solar wind during the time t of the 
duration of a fiight, it is worthwhile to consider y — s/t as a measure of the physical 
variable within a fiight. From the same data, one can obtain distributions presented 
in Fig. 5. By choosing 




(42) 
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where a2 can be computed in the cases of the solar wind and the auroral electro jet 
index, and this yields 

asw ~ 2.0 , oal ~ 2.1 . (43) 
Then the results (7) and (43) give 

4'^i/asw ~ P^l/aAL ~ 0.45 - 0.5 
4w/«sw~ 0.75 -1.1 

Al/«al ~ 0.62 - 0.86 (44) 

The results (44) show that processes SW and AL are very close to each other in zone 
1, related to small scale flights, and fairly different in zone 3 that is related to large 
scale flights. 

In addition, for short flights the corresponding value of /i is close to the normal 
one and even could be related to subdiffusion (// < 1) while for the long flights // > 1 
and is related to superdiffusion processes if the variable y — s/t is the right one with 
respect to the Kolmogorov condition. 



7. NUMERICAL STUDIES OF FRACTIONAL KINETIC EQUATIONS 

In this section numerical results from the investigation of two classes of kinetic 
equations discussed earlier are presented. The flrst is the fractional kinetic equation 
(FKE) as defined earlier, Eq. (33), with (3 = 1 and a = 1.6. The second is the 
Fokker-Planck type equation, Eq. (18), but with a specific power- law dependence for 
the diffusion coefficient. For the sake of definiteness, the two equations are 

dF{y,t) ^ d^F{y,t) 

dt d\y\'^ ^ ' 

dt 9?/^'^' dy ^ ^ ^ 

For these studies fi = 2 — ato make the mean-squared displacement scale similarly 

with time for the two equations. For the FKE, a specific choice of the fractional 
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Figure 5: Distributions P{s/t) for SW (a) and AL (c) and their tails (b) and (d). 
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Figure 6: logio{F) versus logiQ(y) for t — 0.5 (solid), t — 1.0 (dot), t — 1.5 (dash), 
asymptotic theory (dot-dash). 

derivative in y has been made. For this choice, the Fourier transform representation 
for the dependent variable yields 



For Eq. (46) an implicit finite difference scheme yields a tridiagonal equation which 
can be readily inverted to obtain the distribution function F{y, t) at different in- 
stants of time. Although the scheme is a stable one, the time-stepping restriction is 
determined by the requirement on accuracy. 

Fig. 6 is a log-log plot for the solutions of Eq. (45) for a— 1.6 at different instants 
of time. The different lines are for three instants of time namely, t — 0.5 (solid), 
t = 1.0 (dot), t = 1.5 (dash). The last line (dot-dash) is the theoretically determined 
asymptotic solution derived in Sec. 4 (Eq. 40). 




(47) 
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Figure 7: logio(-F) versus logio(y) for t = 2.0 (solid), t = 6.0 (dot), t = 10.0 (dash). 

In Fig. 7 is shown the log-log plot of the solutions of Eq. 46 for /i—OA also at 
different instants of time. Once again the different hues are ai t — 2 (sohd), t — 6 
(dot) and t = 10 (dash). For this case the distribution does not display any power- law 
type behavior. 

Furthermore in Fig. 8, the normalized even moments (to 2n — 10) are plotted for 
the Gaussian distribution (diamonds), Eq (46) (at t — 8) (squares) and for the FKE 
ai t — 1.5 (circles). Here again it is clear that for Eq. (46) the moments are finite 
and comparable to the Gaussian case. However for the FKE the moments are large 
and dependent on the size of the computation box (i.e. maximum value of y used in 
the computation) and in principle are infinite as expected, based on the arguments 
presented in Sec. 4. 

8. CONCLUSION 

In the solar wind-magnetosphere interaction, multiscale features coexist along with 
the global or coherent features, and have been studied extensively using nonlinear dy- 
namical techniques. The detailed properties of these phenomena are studied using 
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Figure 8: < |/^" > / < >" versus n for Gaussian distribution function( diamonds), 
distribution function obtained from Fokker-Planck like equation (Eq. (46)) (squares) 
and distribution function obtained from the fractional kinetic equation (Eq. (45)) 
(triangles) . 
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Levy flights and fractional kinetics models. In the solar wind - magnetosphere cou- 
pling, a technique to analyze the Levy-type processes is applied to the time series data 
of the solar wind electric field and the auroral electrojct index. The probability dis- 
tribution function of the flights show similarities and differences, and provides a new 
insight into the origin of the multiscale behavior through the different values of the 
scaling indices. In a complementary approach, the fractional kinetic equations, which 
uses fractional derivatives to represent the complexity, provide a suitable mathemat- 
ical framework for the multiscale behavior. Unlike the usual diffusion equations, the 
solutions of these equations yield non-convergent moments, showing its multiscale 
features. The origin of multiscale behavior is studied by using a fractional kinetic 
equation and a diffusion equation with a space dependent diffusion coefficient. Nu- 
merical solutions of these equations are used to analyze the nature of the equations 
by computing the moments. The fractional kinetic equations yield solutions with 
large moments and are divergent. On the other hand the solutions of the diffusion 
equation have convergent moments similar to those of Gaussian distributions. These 
results lead to the conclusion that the fractional kinetic equations are suitable models 
of multiscale phenomena. 

In the the study of the solar wind by Hnat et al [5], all the probability distribu- 
tion functions, e. g., their Figs. 3, 4, and 5, have infinite first moment because of the 
small Hurst exponent. However this does not follow from their equations (Eqs.(2) and 
(3)). This inconsistency raises the question of what leads to multiscale properties. 
In the studies presented in this paper, only the fractional kinetic equations lead to 
the formation of tails in the distribution functions, and thus can be considered as 
appropriate models of multiscale processes. It should however be noted that the pa- 
rameters defining fractional kinetic equation could be improved as more data become 
available and with a better understanding of the magnetosphere dynamics. 
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